from problembase import ProblemBase
from dolfin import RectangleMesh, Expression, plot
from math import sqrt,pi,sin,cos,asin
from mikael.Zalesak_mesh import get_phi
from mikael.LagrangianParticles import zalesak

class Zalesak(ProblemBase):
    'Zalesak disk.'
    def __init__(self,options):
        ProblemBase.__init__(self,options)

        self.center = (0.5,0.75)   # center of the disk
        self.r  = 0.15             # disk's radius
        self.w = 0.05              # width of the notch
        self.l = 0.25              # length of the notch

        self.LL = (0,0)
        self.UR = (1,1)
        self.domain = RectangleMesh(self.LL[0],self.LL[1],self.UR[0],self.UR[1],2*self.N,2*self.N)
        
        self.a = -pi
        self.b = 0.5
        self.vStrings = ('-1*%g*(x[1]-%g)' % (self.a,self.b),'%g*(x[0]-%g)' % (self.a,self.b))     # string for velocity compononets    
        self.vIsTimeDependent = False               
        self.hasExactSolution = True
        
        self.N = 50 # number of particles on the disk 
        self.points = zalesak(self.center,radius=self.r,width=self.w,slot_length=self.l,N=self.N) # markers seeded on the disk 
       
        self.phi_ = get_phi(self,0)  # call mikael's function to get the initial solution   
        self.T = 5      
        self.activeCols = [0,1,2,3,4,5,6,7]
        ProblemBase.register_variables(self)
        self.ibcExpr =  self.inflow_boundary_value  # use signed distance from the disk
                        #lambda t : get_phi(self,t) # use constructed exact solution

       
        theta = asin(self.w/2./self.r)
        self.L = self.w + 2*self.l + 2*self.r*(pi-theta)

    def exact_volume(self):
        '''Return problem specific exact volume.'''
        w = self.w
        l = self.l
        r = self.r
        return pi*r**2 - l*w - r**2*asin(w/2./r) + sqrt(r**2-w**2/4.)*w/2.
       
    def exact_solution(self,t):
        '''Signed distance function computed using markers interpolated to finer mesh.'''
        phi_exact = get_phi(self,t)  # on 100x100 mesh, #TODO make higher 
        return phi_exact

    def inflow_boundary_value(self,t):
        '''Signed distance function from the full disk.'''
        a = self.a
        b = self.b
        cx = self.center[0]
        cy = self.center[1]
        return Expression('sqrt(\
                             (b+cos(a*t)*x[0]+sin(a*t)*x[1]-b*cos(a*t)-b*sin(a*t)-cx)*(b+cos(a*t)*x[0]+sin(a*t)*x[1]-b*cos(a*t)-b*sin(a*t)-cx)+\
                             (b-sin(a*t)*x[0]+cos(a*t)*x[1]-b*cos(a*t)+b*sin(a*t)-cy)*(b-sin(a*t)*x[0]+cos(a*t)*x[1]-b*cos(a*t)+b*sin(a*t)-cy)\
                               )-r',r=self.r,cx=cx,cy=cy,a=a,b=b,t=t\
                         )
